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SUMMARY 


Assay  data  are  often  fitted  by  a  nonlinear  regression  model 
incorporating  heterogeneity  of  variance.  Typically,  the  standard  deviation 
of  the  response  is  taken  to  be  proportional  to  a  power  0  of  the  mean.  There 
is  considerable  empirical  evidence  suggesting  that  for  assays  of  a 
reasonable  size,  how  one  estimates  the  parameter  8  does  not  greatly  affect 
how  well  one  estimates  the  mean  regression  function.  An  additional 
component  of  assay  analysis  is  the  estimation  of  auxiliary  constructs  such 
as  the  minimum  detectable  concentration,  for  which  many  definitions  exist; 
we  focus  on  one  such  definition.  The  minimum  detectable  concentration 
depends  both  on  0  and  the  mean  regression  function.  We  compare  standard 
methods  of  estimating  the  parameter  0  due  to  Rodbard  (1978).  Raab  (1981a) 
Sadler  &  Smith  (1985)  and  Carroll  and  Ruppert  (1982b).  When  duplicate 
counts  are  taken  at  each  concentration,  the  first  method  is  only  20% 
efficient  asymptotically  in  comparison  to  the  fourth  for  normal  data,  and  in 
an  example  the  resulting  estimate  of  the  minimum  detectable  concentration  is 
asymptotically  3.7  times  more  variable.  Less  dramatic  results  obtain  for 
the  second  and  third  estimators  compared  to  the  fourth.  Simulation  results 
and  an  example  support  the  asymptotic  theory.  The  results  have  implications 
in  applications  other  than  the  assay  problem  in  which  heterogeneity  of 
variance  and  issues  of  calibration  arise. 

I 

Some  key  words:  Calibration;  Generalized  Least  Squares;  Heteroscedasticity; 

Prediction;  Weighted  Least  Squares. 


1.  INTRODUCTION 


Recent  work  In  the  analysis  of  assay  data  In  the  clinical  and  biological 
sciences  suggests  that  .hese  data  can  be  markedly  heteroscedastic.  In 
radioimmunoassay,  this  characteristic  has  been  observed  repeatedly  and 
incorporated  into  the  analysis  as  discussed  by  Finney  (1976),  Rodbard  (1978), 
Tiede  &  Pagano  (1979),  Raab  (1981a,b),  Butt  (1984)  apd  Sadler  &  Smith  (1985). 
Such  analyses  are  for  the  most  part  special  cases  of  the  heteroscedastic 
nonlinear  regression  model.  Specifically,  we  observe  independent  count:.  ^  at 
concentrations  Xj  for  i  =  1.....N  and  j  =  l,..,Mj  with  means  and  variances 
given  by 

E  (Y^)  =  =  f(xrP)  ;  var^)  =  (ag^./J.e)}2.  (1.1) 

where  p  is  the  unknown  regression  parameter  vector  of  length  p,  g  is  the 
variance  function,  and  0  is  the  structural  variance  parameter.  A  standard 
model  for  the  mean  in  a  radioimmunoassay  is  the  four  parameter  log-logistic 
model 


f(x.P)  =  Pj  +  ( P2  ~  Pj)/[ 1  +  exp(P4(log  x  -  P3 )}].  (1.2) 

Almost  without  exception,  the  variances  have  been  modeled  as  functions  of  the 
mean  response,  usually  either  as  a  quadratic  or  as  a  power  of  the  mean,  e.g., 

Oj  =  Standard  deviation  of  Y^j  =  crgfXj.p.B)  =  affXj.p )®.  (1-3) 

j 

i 

The  fundamental  contribution  of  Rodbard  and  other  workers  has  been  to 
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incorporate  the  heterogeneity  into  the  analysis,  thus  improving  the  quality  of 
statistical  analysis. 

Fitting  model  (1.1)  in  assays  is  discussed  by  many  authors,  see  Rodbard  & 
Frazier  (1975),  Raab  (1981a),  and  Finney  (1976).  The  most  common  method  of 
estimating  0  is  generalized  least  squares,  in  which  one  forms  estimates  of  the 
a j  and  then  estimates  0  by  weighted  least  squares.  For  estimating  the  weights, 
both  formal  estimation  procedures  for  estimating  0  (Rodbard. &  Frazier,  1975  and 
Raab,  1981a)  and  the  method  of  setting  6  to  a  predetermined  fixed  value  based 
on  experience  with  a  particular  assay  (Finney,  1976)  have  been  advocated.  We 
consider  estimation  of  the  variance  parameter  to  be  an  important  problem 
independently  of  whether  or  not  0  is  set  to  a  fixed  value  in  curve  fitting  and 
subsequent  analyses.  For  example,  estimates  of  0  from  developmental  and 
validation  runs  of  an  assay  are  commonly  used  to  establish  a  fixed  value,  or 
range  of  acceptable  values,  for  0.  as  in  the  example  of  Section  5.  Also, 
routine  calculation  of  0  after  assay  implementation  provides  information  that 
can  be  helpful  in  monitoring  for  assay  changes.  In  the  development, 
production,  and  quality  control  of  a  pharmaceutical  product,  for  example, 
assays  characteristically  Include  fewer  test  samples  and  greater  numbers  of 
known  concentrations  than  in  many  clinical  applications.  This  is  particularly 
true  when  the  substance  of  interest  is  at  very  low  levels  in  the  final  product. 
In  this  case,  it  is  practical  to  estimate  0,  as  well  as  limits  of  assay 
reliability  as  discussed  below.  Extreme  deviations  from  historical  values  of 
these  parameters  may  be  considered  when  evaluating  the  acceptability  of  an 
assay.  The  intention  of  this  paper  is  to  demonstrate  that  how  well  one 
estimates  the  variance  function  can  be  crucial  in  determining  the  properties  of 
analyses  based  on  the  fitted  curve.  Qualitatively.  then,  how  0  is 
characterized,  whether  by  formal  estimation  or  as  fixed,  is  an  important  issue. 
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When  8  is  formally  estimated,  as  discussed  by  Jobson  &  Fuller  (1980)  and 
Carroll  &  Ruppert  (1982b),  j3  has  the  same  asymptotic  normal  distribution  as  the 
weighted  least  squares  estimator  with  known  a The  first  order  asymptotics 
can  be  optimistic,  but  in  our  experience  for  RIA  and  ELISA  assays,  the 
asymptotics  are  often  reasonable.  Thus,  if  our  only  interest  is  to  estimate  P, 
in  many  assays  the  method  of  estimating  the  variance  function  may  not  be 
crucial.  However,  the  assay  problem  does  not  always  stop,  with  estimating  P. 
but  also  addresses  issues  of  calibration,  such  as  confidence  intervals  for  a 
true  given  a  new  Y  ,  the  classic  calibration  problem,  and  determining  the 
sensitivity  of  the  assay  using  such  concepts  as  the  minimum  detectable 
concentration  of  Rodbard  (1978)  and  the  critical  level,  detection  level  and 
determination  limit  of  Oppenheimer,  et  al.  (1983).  With  the  increase  in 
laboratory  automation  and  the  use  of  computers,  routine  estimation  of  assay 
detection  limits  is  becoming  commonplace.  As  noted  in  practice  by  Oppenheimer, 
et  al.  (1983)  and  others,  a  feature  of  these  calibration  problems  is  that 
efficiency  of  estimation  is  essentially  determined  by  how  well  one  estimates  0. 
This  paper  provides  theoretical  and  empirical  evidence  to  justify  this  claim. 
The  qualitative  implications  of  this  result  are  of  interest  regardless  of  how 
estimation  of  the  variance  structure  is  performed.  For  example,  Finney  (1978, 
p.  342)  describes  how  different  assumptions  about  the  value  of  0  can  markedly 
affect  the  length  of  confidence  intervals  in  potency  estimation  for  RIA  assays. 

To  the  best'  of  our  knowledge,  our  paper  is  one  of  the  first  which  shows 
explicitly  that  how  one  estimates  the  structural  variance  parameter  0  can  be 
Important  in  determining  the  behavior  of  estimates  of  interesting  quantities; 
this  conclusion  extends  beyond  the  realm  of  the  assay  problem.  Far  from  being 
only  a  nuisance  parameter,  0  has  an  important  role  in  the  analysis  of 
calibration  and  prediction  problems.  The  implications  are  of  interest  when  0 


is  important  in  itself,  as  in.  for  example,  off  line  quality  control,  see  Box 
(1987)  who  discusses  finding  the  levels  of  x  which  give  minimum  variance 
subject  to  a  constraint  on  the  mean. 


For  our  development,  instead  of  pursuing  a  fully  general  theory,  we  focus 
on  the  determination  of  minimum  detectable  concentration.  There  is  no  unique 
definition  of  this  concept,  and  for  illustration  we  pick  one  of  the  possible 
candidates  (Rodbard,  1978). 

Definition.  Let  Y(x,M)  be  the  mean  response  based  on  M  replicates  at 
concentration  level  x.  taken  independently  of  the  calibration  data  set  (Y.  } 

*  J 

and  let  f (0,0)  be  the  expected  response  at  zero  concentration.  The  minimum 

detectable  concentration  x  at  level  (1-a)  is  the  smallest  concentration  x  for 

c 

which 


pr(Y(x,M)  l  f(0 .0)}  >  1  -  a. 


(1.4) 


Other  definitions  replace  0  with  0  in  (1.4).  If  t(a.N-p)  is  the  (l-a)th 
percentile  of  the  t-distribution  with  N-p  degrees  of  freedom,  the  usual 

A 

estimate  x  of  x  satisfies 
c  c 


(f(xci)  -  f(O,0)}2  =  {t(a.N-p)}2{a2g2(xci.0)/M  +  v[f(O,0)]}. 


(1.5) 


^  A  AA 

where  v[f(O,0)]  is  tin  estimate  of  the  variance  of  f(0 ,0)  and  a  is  the  usual 
mean  squared  error  from  the  weighted  fit.  To  be  precise  one  would  replace  the 
t-percentage  point  with  a  correction  based  on  the  limit  distribution  of  the 
estimates  of  (0,0, o),  but  this  has  not  been  followed  in  practice  since  the 
limit  distribution  has  been  unknown  and  the  effect  is  asymptotically 
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unimportant.  A  heuristic  discussion  of  how  estimation  of  the  variance 
structure  plays  a  role  in  the  determination  of  calibration  quantities  can  be 
found  in  Davidian  &.  Carroll  (1987). 

In  Section  2.  we  discuss  three  standard  methods  for  estimating  the 
variance  parameter  0.  In  Section  3  we  show  that,  under  relatively  general 
conditions,  two  of  these  can  be  quite  a  bit  less  efficient  than  the  third  and 
discuss  how  this  difference  translates  theoretically  to  the  minimum  detectable 
concentration  problem.  In  assays,  along  with  standard  observations  at  known 
concentrations,  there  are  also  test  samples  composed  of  additional  observations 
at  unknown  concentrations.  The  first  two  estimators  can  incorporate 
information  from  both  standards  and  unknowns  while  the  third  is  equipped  only 
to  use  Information  from  standards.  We  discuss  how  to  combine  estimators  to 
accomodate  both  types  of  data  with  an  increase  in  efficiency.  In  Sections  4 
and  5.  we  discuss  a  small  Monte-Carlo  study  and  an  example  to  illustrate  the 
results.  The  key  conclusion  is  that  how  one  estimates  0  can  affect  the 
relative  efficiency  of  estimated  quantities  useful  in  the  calibration  of 
assays . 


2.  METHODS  OF  ESTIMATING  MEAN  AND  VARIANCE  PARAMTETERS 
The  problem  of  estimating  0  in  models  (1.1)  and  (1.3)  has  been  discussed 
in  many  places  in  the  literature.  See,  for  example,  Rodbard  (1978),  Jobson  & 
Fuller  (1980),  RAab  (1981a),  Carroll  &  Ruppert  (1982b)  and  for  a  general  theory 
Davidian  &  Carroll  (1987).  We  focus  our  attention  on  three  methods,  two  of 
which  were  proposed  in  the  assay  context.  For  simplicity,  we  discuss  only  the 
case  of  equal  replication  Mj  =  M  j>  2.  The  first  two  methods  require  some 
replication  and  do  not  use  the  form  of  the  mean  response  so  that  standards  and 
unknowns  both  may  be  considered. 
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2.1  Loe- linear! zed  estimation 


Model  (1.3)  implies  upon  taking  logarithms  that  the  log  standard  deviation 

_  2 

is  linear  in  the  log  mean  with  slope  0.  If  (Yj.Sj  )  are  the  within 

concentration  sample  means  and  variances,  this  suggests  estimation  of  0  as  the 

- 

slope  from  regressing  log  Sj  on  log  Yj.  If  we  denote  this  estimate  as  0^. 
Rodbard  (1978)  suggests  forming  estimated  standard  deviations  as  the  sample 
means  to  the  power  0;  T  .  and  then  applying  weighted  least  squares  to  estimate  p. 

2.2  Modified  maximum  likelihood 

Assuming  independence,  Raab  (1981a)  suggests  estimating  0  without  making 
any  assumptions  about  the  form  of  the  mean  function.  Raab  proposes  estimation 
of  0  by  joint  maximization  of  the  ''modified"  normal  likelihood 

N  M 

U  {2mcr2g(p i.0))(M'l)/2eXp[-  I  (Y  .  .-p  )2/(2o2g2(p .0)}]  (2.1) 

i=l  1  J=1  U  1  1 

in  the  parameters  a,  0,  p^ . p^,  where  we  have  written  g (x^.p.d)  as  g(pit0) 

to  emphasize  the  dependence  of  the  variance  function  on  the  mean  response. 
Estimation  of  P  may  now  proceed  via  weighted  least  squares  in  a  fashion 
analogous  to  the  log-linearized  method.  Sadler  &  Smith  (1985)  maximize  (2.1) 
in  a  and  0  but  with  Pj  estimated  by  Yj .  Their  estimate  of  0  is  easier  to 
compute  and  asymptotically  equivalent  to  that  of  Raab  under  the  conditions  of 
Section  3. 

1 

2.3  Pseudo-1 ikel lhood 

For  given  p  the  pseudo-1 ike lihood  estimator  of  0  is  the  normal  theory 
maximum  likelihood  estimate,  maximizing 

-M2  log{g(x  p,0))  -  (N/2) log[N  2  2  (Y  -f  (x.  .P)}2g(x.  ./3,0)"2]. 

i=l  i=l  j=l  1  1 
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see  Carroll  &  Ruppert  (1982b).  To  estimate  0  and  P  jointly  one  can  (1)  set  P  = 
unweighted  least  squares;  (ii)  estimate  0  by  pseudo-likelihood;  (111)  form 

2  A  A 

estimated  variances  g  (Xj.0.0);  (iv)  re-estlmate  P  by  weighted  least  squares; 
(v)  iterate  (li)  -  (iv)  one  or  more  times.  The  number  of  cycles  *0  of  this 
algorithm  is  the  number  of  times  one  performs  step  (iv).  One  can  do  step  (ii) 
by  direct  maximization  or  by  weighted  least  squares  as  in.Davidian  &  Carroll 
(1987). 

Pseudo- likelihood  requires  no  replication  and  easily  copes  with  unequal 
replication.  Since  pseudo- likelihood  depends  on  the  form  of  the  mean  function, 
only  data  from  standards  may  be  used  to  construct  the  estimate. 

2.4  Other  methods 

Other  methods  have  been  proposed;  see  Jobson  &  Fuller  (1980)  and  Box  &. 
Hill  (1974).  Robust  variance  function  estimation  methods  have  also  been 
developed  by  Carroll  &  Ruppert  (1982b)  and  Giltinan,  Carroll  &  Ruppert  (1986). 

A  final  method  of  jointly  estimating  (/?.©)  is  normal  theory  maximum 
likelihood.  There  are  important  Issues  of  robustness  which  complicate  routine 
use  of  this  method,  see  McCullagh  (1983),  Carroll  &  Ruppert  (1982a)  and 
Davidian  «  Carroll  (1987)  for  further  discussion.  For  assay  data, 
pseudo-likelihood  and  maximum  likelihood  estimates  of  0  have  similar  asymptotic 
behavior,  and  we  use  the  former  largely  for  its  ease  of  calculation. 

At  this  poiht  we  address  briefly  issues  related  to  assay  design.  It  may 
be  argued  that  estimation  of  the  variance  structure  ignores  the  difficulties 
introduced  by  lack  of  randomization,  in  particular  that  adjacent  positioning  of 
replicates  will  likely  lead  to  underestimation  of  the  variance.  This  is  an 
important  but  difficult  issue  which  involves  problems  of  variance  component 
estimation  and  correlated  observations,  and  we  consider  it  beyond  the  scope  of 
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the  discussion  here.  It  is  worthwhile  to  note  that  the  pseudo- 1  ike  1  ihood 
estimator  includes  lack-of-fit  error  as  well  as  replication  error  while  the 
log- linear! zed  estimator  does  not  Include  the  former.  The  two  methods  do  not 
estimate  the  same  quantity  without  the  proper  assumptions. 


3.  ASYMPTOTIC  THEORY 

A 

The  asymptotic  theory  of  the  log-linearized  estimator.  0^,,  is  complicated 
because  regressing  log  Si  on  log  Yj  is  not  a  standard  linear  regression 
problem.  Likewise,  the  asymptotic  theory  for  the  modified  maximum  likelihood 
estimator  *s  complicated  because  the  dimension  of  the  parameter  space 

increases  with  N.  Both  problems  may  be  treated  as  nonlinear 


errors-in-variables  problems  as  in  Wolter  &  Fuller  (1982)  and  Stefanski  & 

__  A 

Carroll  (1985).  The  error  in  estimating  p 4  by  Yj  in  8t  r  or  by  the  joint 

A  A 

estimator  p^  in  causes  these  estimators  to  be  biased  asymptotically.  This 

bias  is  typically  negligible,  because  in  most  assays  the  parameter  a  in  (1.1) 
is  quite  small.  It  thus  makes  sense  to  define  an  asymptotic  theory  where  the 


sample  size  Ng  =  NM  becomes  large  and  a  simultaneously  is  small.  Because  in 
most  assays  the  number  of  replicates  M  is  small,  we  let  N  -*  «  and  a  -*  0  while 


keeping  M  fixed;  Raab  (1981a)  suggests  that  M  =■  2  is  most  common.  The  approach 
of  letting  N  -*  •  and  a  -*  0  simultaneously  is  dictated  by  the  problems  of 

A  A  /v 

studying  0LL  and  SyjyjL-  The  pseudo-likelihood  estimator  0p^  has  routine 
asymptotic  theor^  even  for  fixed  a. 


The  asymptotic  distributions  of  these  estimators  for  0  can  be  obtained 


from  the  general  theory  of  Davidian  &  Carroll  (1987).  Throughout  this 
discussion  we  focus  on  model  (1.3).  Define  =  {Yjj  -  f (Xj ,P))/(of (xA ,p)0} . 
Vj  =  log  ffXj.p).  q2  =  (M-l)“1I(elj-ri)2,  and  a2  =  lim  (N-l )_1I(vi-7)2. 
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THEOREM  1.  As  N  -*  «  and  a  -*  0  simultaneously  and  N  a  =  0(1),  if  the  random 
variables  (a  }  are  symmetric  and  independent  and  identically  distributed,  then 

*  J 

A  A  A* 

®LL*  ®MML  ®PL  are  normally  distributed  with  mean  0  and 

variances  var{log  q1)/(-4Na^).  var(q^)/(4Na^) ,  and  var(e^j)/(4NMo^) . 

respectively. 


Here,  the  symmetry  condition  is  necessary  to  ensure  that  the  asymptotic 

A  A 

distributions  of  0T  t  and  zero  mean;  symmetry  is  unnecessary  for  the 

result  for  0p^.  Since  var(qi)  =  {var(e^j)/M}  +  [2/{M(M-l)}] ,  has 

uniformly  larger  asymptotic  variance  than  0p^  for  all  M  £  2  regardless  of  the 
distribution  of  the  {e^}.  for  normal  data,  the  efficiency  of  modified  maximum 
likelihood  relative  to  pseudo-likelihood  is  (M-l)/M.  A  simple  calculation 
shows  in  this  case  that  the  efficiencies  of  the  log-linearized  method  relative 
to  modified  maximum  likelihood  are  40. 5%,  60.8%,  71.3%,  and  89.3%  for  M  =  2.  3, 

A  A 

4.  and  10.  respectively,  so  that  the  efficiency  of  8,  T  relative  to  0p^  is  only 
20.3%  for  M  =  2  and  53.5%  for  M  =  4.  In  our  experience,  these  numbers  slightly 
exaggerate  the  inefficiency  of  8| t  relative  to  pseudo- likelihood,  especially  in 
assay  problems  where  N  is  rather  small,  as  in  the  Monte  Carlo  study  of  Section 

A  A 

4.  The  asymptotic  relative  efficiencies  of  8,  T  relative  to  aSree  quite 

well  with  the  efficiencies  of  39%,  62%  and  77%  and  39%,  64%  and  74%  for  M  =  2, 
3  and  4  reported  by  Raab  (1981a)  and  Sadler  &  Smith  (1985),  respectively,  in 
two  Monte-Carlo  ‘studies  with  larger  N.  Davidian  &  Carroll  (1987)  show  that 
under  the  condition  a  -»  0  and  (1.1).  pseudo- likelihood  and  normal  theory 
maximum  likelihood  are  asymptotically  equivalent.  Thus,  while  the  inefficiency 
of  the  log- linearized  and  modified  maximum  likelihood  methods  is  not  surprising 
for  normal  data,  Theorem  1  shows  the  large  extent  to  which  the  log-1 inear i zed 
method  can  be  inefficient  for  small  M  under  normality  as  well  as  the  fact  that 


modified  maximum  likelihood  is  inefficient  for  all  distributions. 

A  benefit  of  the  log- linear! zed  and  modified  likelihood  methods  is  that 
they  can  incorporate  the  responses  at  unknown  concentrations  to  estimate  0 
while  pseudo- likelihood  can  not.  One  can  improve  upon  pseudo- likelihood  to 
take  into  account  unknowns  by  considering  a  weighted  estimator  of  0  which 

-  A 

combines  0p^  and  one  of  the  others  in  the  obvious  way.  The  increase  in 
efficiency  of  such  an  estimator  over  the  others  alone  will  depend  on  the  amount 
of  information  available  from  test  samples  relative  to  standard  samples.  In 
this  discussion  we  have  restricted  our  investigation  to  standard  observations 
because  the  improvement  in  efficiency  by  pseudo-likelihood  comes  from  these 
data. 

For  the  minimum  detectable  concentration,  note  that  in  (1.5)  the  term 

-i 

v{f(0 ./?)}  is  of  the  order  (NM)  and  is  hence  small  relative  to  all  the  other 
terms.  Of  course,  for  normally  distributed  data,  the  solution  to  (1.5)  is  the 
quantity  x*.  where 

0  =  (z(a)}2a2f20(x*.p)/M  -  (f(x*  P)  -  f(0.p)}2.  (3.1) 

and  z(a)  is  the  (1  -  a)**1  percentile  point  of  the  standard  normal  distribution. 

Here  is  the  major  result,  the  technical  details  for  which  are  given  in  the 
appendix.  Define  dQ  =  log  f(O.P)  -  lim  N_12  log  ((x^.p). 

i 

A  A  A 

THEOREM  2.  Let  x  (LL),  x  (MML)  and  x  (PL)  denote  the  estimated  minimum 

c  c  c 

detectable  concentrations  using  the  log- linearized  estimate  0^,  the  modified 

A 

maximum  likelihood  estimate  the  pseudo- likelihood  estimate  0p^, 

respectively.  Then  under  regularity  conditions  and  the  conditions  of  Theorem 

A  A 

1.  there  is  a  constant  AQ  and  a  sequence  b^  for  which  x  (LL).  x  (MML)  and 
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xc(PL)  are  asymptotically  normal  with  mean  xq  and  variances  Vt ;  cN ,  Vjt^c^  and 
VPLCN*  resPectively*  where  =  a2/(Nb^A2).  =  var(e2j)+var(log  q2)d^M/o2) . 

VMLL  "  var (£^ j )+var (q2)d^M/o2 ,  and  VpL  =  varfe^  ){l+d2/a2}  . 

From  the  above,  the  ordering  in  efficiency  of  estimated  minimum  detectable 

concentration  is  the  same  as  the  ordering  of  estimators  of  6.  The  asymptotic 

relative  efficiency  of  the  modified  maximum  likelihood  estimate  of  minimum 

detectable  concentration  relative  to  the  pseudo- likelihood  estimate  is 

var(fe2j)(a2+dQ)/{var(e2j)(a2+d2)+2d2/(M-l)} .  which  is  less  than  1  for  all  M 

2 

regardless  of  the  value  of  var(e  .).  The  asymptotic  relative  efficiency  of  the 

*  J 

log- linearized  estimate  of  minimum  detectable  concentration  to  the 

2  2  2  2  2  2  2 

pseudo-likelihood  estimate  is  var(£ij)(ov+dQ)/{ovvar(eij)+Hd0var(log  qj)}. 
Thus,  pseudo- likelihood  is  favored  over  modified  maximum  likelihood  for  any 
underlying  distribution  and  is  favored  over  the  log- linear! zed  method  at  least 
when  the  data  are  approximately  normal.  For  distributions  other  than  normal, 
calculations  with  other  symmetric  distributions  such  as  double  exponential  and 

A 

various  contaminated  normal  distributions  show  very  few  cases  where  6{  t  is  more 

A 

efficient  than  0p^.  The  numerical  efficiencies  depend  on  the  logarithm  of  the 

2  2 

true  means  through  d^  and  a ■  .  For  example,  in  the  simulation  discussed  in  the 
next  section,  the  asymptotic  relative  efficiency  of  the  log- linear! zed  estimate 
is  26.6%  for  M  =  2  and  62. IX  for  M  =  4. 

The  asymptotic  theory  thus  confirms  that  inefficiencies  in  estimating  the 
variance  parameter  0  translate  into  inefficiencies  for  estimating  the  minimum 
detectable  concentration. 


4.  A  SIMULATION 


To  check  the  qualitative  nature  of  the  asymptotic  theory,  we  ran  a  small 
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simulation  based  on  an  ELISA  assay.  The  responses  Y^j  were  normally 
distributed  with  mean  and  variance  satisfying  (1.2),  (1.3),  where  0^  =  30.5754, 
02  =  1.9173,  03  =  1.7235,  04  =  0.9730.  0  =  0.7  an d  a  =  0.0848.  The  12 
concentrations  were  chosen  from  the  example  of  Section  5  to  represent  the 
number  of  levels  typical  for  assay  data.  These  were  every  other  concentration 
beginning  with  0.000  in  Table  1.  For  each  of  500  data  sets  for  M  =  2  and  4,  0 
was  estimated  by  unweighted  least  squares  and  by  generalized  least  squares  with 
0  estimated  by  the  log- 11 near! zed,  pseudo-likelihood,  and  modified  maximum 
likelihood  methods,  with  weighting  as  in  Section  2.  We  used  the  estimate  of 
Sadler  &  Smith  (1985)  described  in  Section  2  in  place  of  modified  maximum 
likelihood  and  computed  the  pseudo- likelihood  estimate  for  *6=2  cycles  of  the 
algorithm  as  suggested  by  the  results  of  Davldian  &  Carroll  (1987).  The 
estimates  of  0  were  constrained  to  lie  in  the  interval  0  £  0  £  1.50. 

For  M  =  2.  the  Monte  Carlo  biases  and  standard  errors,  in  parentheses,  for 
0  were  -0.0371  (0.2848)  for  the  log- linearized  method,  -0.0162  (0.2146)  for 
modified  maximum  likelihood  and  -0.0107  (0.1799)  for  pseudo-likelihood.  For  M 
=  4  these  were  -0.0026  (0.1504),  -0.0009  (0.1316)  and  -0.0012  (0.1247), 
respectively.  Thus,  Monte  Carlo  efficiencies  relative  to  pseudo- likelihood 
based  on  Monte  Carlo  mean  squared  errors  were  for  M  =  2  (M  =  4)  39.4%  (68.7%) 
for  the  log- linearized  method  and  70.1%  (89.8%)  for  modified  maximum 
likelihood,  compared  to  theoretical  values  of  20.3%  (53.5%)  and  50.0%  (75.0%), 
respectively.  The  asymptotics  tend,  to  exaggerate  the  loss  of  efficiency; 
nonetheless,  this  example  indicates  that  the  pseudo- likelihood  estimator  can  be 
in  some  circumstances  a  considerable  Improvement  over  the  other  methods. 

For  the  minimum  detectable  concentration  we  chose  a  =  0.05.  For  all 
methods,  a  Monte  Carlo  measure  of  (1.4)  based  on  generating  means  of  M  new 

A 

responses  at  xc  for  each  data  set  showed  that  this  requirement  was  satisfied; 
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rather  than  0.95,  every  case  was  at  least  0.97.  Mean  values  of  estimated 
minimum  detectable  concentrations  multiplied  by  100  and  their  Monte  Carlo 
variances  relative  to  pseudo- like  11  hood  (in  parentheses)  are  reported  below'. 
Monte  Carlo  Results  for  Minimum  Detectable  Concentration 


M  =  2 


M  =  4 


Least  Squares 
Log- 1 inear 1 zed 


Modified  Maximum 


14.023  (10. 52) 
5.102  (33.22) 


9.747  (8.82) 
3.033  (74.92) 


Likelihood 


4.685  (52.82) 


Pseudo- likelihood  4.313  (100.02) 


3.010  (84.62) 
2.928  (100.02) 


The  relatively  poor  behavior  of  unweighted  least  squares  is  evident.  To  quote 
from  Oppenheimer,  et  al.  (1983)5  "Rather  dramatic  differences  have  been 
observed  depending  on  whether  a  valid  weighted  or  inappropriate  unweighted 
analysis  is  used."  The  mean  minimum  detectable  concentration  for  the 
log- linear! zed  method  was  182  larger  than  for  pseudo- likelihood  method  for  M  = 
2  and  42  larger  for  M  =  4.  Whether  raw  numerical  difference  is  of  any 
practical  consequence  will  depend  on  the  context.  The  method  of  estimating  0 
seems  to  have  important  consequences  for  the  variability  of  the  estimate  of 
minimum  detectable  concentration.  The  asymptotics  suggest  efficiencies  of 
26.62  for  the  log- linearized  method  and  58.82  for  modified  maximum  likelihood 
when  M  =  2  and  62.12  and  81.12,  respectively,  when  M  =  4,  so  overall  the 
asymptotics  exaggerate  somewhat  the  loss  of  efficiency.  We  also  computed  the 
estimates  of  minimum  detectable  concentration  when  0  was  known  and  equal  to 
0.7,  so  that  the  only  difference  between  the  methods  was  whether  the 
generalized  least  squares  weighting  depended  on  sample  means  or  predicted 
values.  In  every  Instance,  the  mean  estimates  of  minimum  detectable 


concentration  and  their  variances  were  very  close,  suggesting  that  the  method 
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used  to  estimate  0  was  responsible  for  the  differences  among  the  estimates  and 
their  variances  reported  above. 

The  qualitative  implication  is  that  the  relative  efficiency  of  the 
estimated  minimum  detectable  concentration  can  be  affected  by  the  algorithm 
used  to  estimate  the  variance  parameter  0. 

5.  AN  EXAMPLE 

We  illustrate  the  implications  of  the  theory  by  an  example  from  the 

development  stage  of  an  ELISA  assay  for  which  additional  concentrations  were 

used  to  characterize  the  mean  and  variance  of  the  standards,  see  Table  1.  The 

analysis  Is  for  illustrative  purposes  only  and  shows  that  the  methods  can  lead 

to  nontrivial ly  different  results  even  with  a  larger  sample  size. 

For  the  full  data  set  and  reduced  data  sets  considering  all  possible 

combinations  of  duplicates  (except  one  set  for  which  an  =  0).  we  computed 

assuming  (1.2)  and  (1.3)  the  pseudo- like 11 hood,  log- 11 near! zed  and  Sadler  & 

Smith  estimates  of  0  and  x  and  that  of  x  based  on  unweighted  least  squares. 

c  c 

The  results  are  given  in  Table  2  and  show  that  the  three  estimates  vary 

greatly.  As  a  crude  measure  of  this  we  computed  the  means  and  standard 

A  A 

deviations  of  0  and  xq  for  the  five  data  sets  obtained  from  duplicates 

(ignoring  the  fact  that  these  are  not  strictly  independent).  Below  we  list 

"relative  efficiencies"  based  on  these  crude  measures: 

"Relative  Efficiencies"  .for  Estimators  of  0  and  x 
- c 

LL  to  PL  MML  to  PL  LL  to  MML 

0  22. 2%  35. IX  63. 2X 

xc  52. 9%  65.9%  80.2% 

Qualitatively,  the  estimators  exhibit  the  behavior  predicted  by  the  theory: 
quantitatively,  the  values  compare  favorably  with  the  theory  given  the  crudity 
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of  the  comparison.  The  example  shows  thai  there  can  be  wide  differences  among 
the  various  methods  for  estimation  of  0  and  minimum  detectable  concentration. 

^  6.  DISCUSSION 

We  have  addressed  the  general  issue  of  estimating  calibration  quantities 
in  assays  which  exhibit  large  amounts  of  heterogeneity.  The  wrong  model 
assumption  of  not  weighting  at  all  leads  to ‘large  decreases  in  the  efficiency 
of  analysis.  Even  when  weighting  is  used,  changes  in  relative  efficiency  occur 
depending  on  the  method  of  estimating  the  variances,  especially  the  parameter  0 
in  (1.3).  The  key  point  is  that  while  for  estimation  of  (5  the  effect  of  how 
one  estimates  the  variance  function  is  only  second  order,  for  estimation  of 
other  quantities  such  as  minimum  detectable  concentration,  the  effect  is  first 
order.  The  implications  of  the  results  are  relevant  in  other  calibration 
contexts  besides  the  assay  problem. 

We  have  had  success  using  the  idea  of  pseudo-likelihood  in  Carroll  fit 
Ruppert  (1982b);  this  method  applies  in  general  heteroscedastic  models  and  is 
easy  to  compute,  although  the  reader  should  be  aware  that  it  is  not  robust 
against  outliers.  In  the  assay  problem,  the  pseudo-likelihood  method  may  be 
combined  with  one  of  the  others  in  the  case  in  which  standard  as  well  as 
unknown  observations  are  present,  resulting  in  an  increase  in  efficiency.  The 
gain  in  efficiency  will  depend  on  the  particular  data  set. 

I 

One  can  also  consider  data  transformation  rather  than  weighting.  The 
transform-bo th-sides  idea  in  Carroll  &  Ruppert  (1984)  applies  to  the  assay 
problem.  Further  results  and  descriptions  of  other  methods  of  variance 
function  estimation  are  given  by  Davldlan  8i  Cbrroll  (1987). 
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APPENDIX.  PROOFS  OF  RESULTS. 

The  analysis  of  the  minimum  detectable  concentration  is  complicated  by  the 

behavior  of  the  derivative  of  f(x,0)  with  respect  to  0  at  x=0,  especially  for 

the  standard  model  (1.2).  Write  f(x,0)  i  h(Tj.0)  .  •  where  tj  =  «(x.0).  vc  = 

«(x*.0)  and  ii  =  «(x  .0).  In  the  model  (1.2).  ti  =  «(x.0)  =  exp(04  log  x) .  We 
c  c  c 

assume  throughout  that  f(O.0)  >  0,  and  that  all  functions  are  sufficiently 
smooth.  Assume  further  that 


<(0,p)  *  0  ;  (A.l) 

d/dn  h(O.0)  =  h^o.0)  /  0  :  (A.2) 

0,0)  =  0  ;  (A. 3) 

If  w  -»  0  and  v  is  a  random  variable  such  that 
p-lim  £(v,0)/£(w,0)  =  1,  then 

p-lim  sup{  |  £^(av+(l-a)w,p)/£^(w,0)  -  1  |  for  0<[a£l  }  =  0  .  (A. 4) 


These  assumptions  are  satisfied  for  the  model  (1.2)  if  04  >  0. 

The  proofs  of  Lemmas  A.2  and  A. 3  are  at  the  end  of  the  appendix.  Let  c  = 

w«)}2. 

1 

1 

if  2 

LEMMA  A.l.  As  ff  ■*  0  for  f(O,0)  >  0.  n  =  oa  +  0(a  )  and  a  = 

C  u  l* 

(c/M)1/2f9(O.0) {d/dn  h(O,0))-1 . 

Proof :  A  Taylor  series  expansion  of  (3.1)  in  n  and  around  zero. 


LEMMA  A. 2.  Assume  that  as  N  -»  ®,  a  -»  0,  (tjc  -  tjc)  =  0p(aN  ).  Define  Aj  = 

2Ma Jd/dr)  h(O./3)}2/(cf20(O./3)}.  Then  as  N  -*  ®.  a  ->  0.  if  N1/2(0  -  6)  =  0  (1). 
c  p 

>  yjl/2^  *w 

A1N  ^c~^c)/a 

1/9  A9  9  9  1/9  A 

=  N 1,Z(a  -  a)/a  +  2{log  f(0./3)}  N1/Z(0  -  6)  +  op(l).  (A. 5) 

LEMMA  A. 3.  Consider  Lemma  A. 2,  Then 

1/9  A9  0  0  -1/9  ^  ^  9  _  i/o  * 

(NM)i/z(a  -  a)/a  =  (NM)  1/Z  2  2  (tf.  -  1)  -  2v(NM )  (0  -  0)  +  o  (1) 

1=1  j=l  J  p 

1/2 

so  that  if  Aq  =  M  Aj  and  d^  is  defined  as  in  Section  3. 

A0Nl/2(^c  "  \)/a  =  (NM)"172^  -  X)  +  2dQ(NM)1/2(0-0)  +  0p(  1 ) . 

PROPOSITION  1  :  The  results  of  Theorem  2  hold  for  tj  ,  where  tj  =  t)  (LL), 

c  c  c 

T)c(MML)  or  nc(PL). 

Proof  of  Proposition  1  :  From  Davldian  &  Carroll  (1987),  using  Theorem  1  we 
have  that 

1/2  A  1  _i/9  ^  o  O  _ 

N  <0LL  "  0)  =  (1/2)  N  2  { loK  qf  -  E(log  qf)}  (v.  -  v)  +  ofl)  : 

1=1  1  11  P 

N 

n1/2(0MML  "  0)  =  <1/2)  N_1/2i21<q?  -  (yi  -  y)  +  Op(l)  ; 

N1/2(0pL  -  0)  =  (1/2)  N"1/2  M"1  2  2  (e  2  -  1)  (v  -  7)  +  o  (1)  . 

1=1  J=1  J  1  p 


so  that  by  Lemmas  A.l  -  A. 3  and  equation  (A. 5),  we  have 
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aoh1/2{Vpl)  "  \)/a 

^  M 

=  (NM)"1/2  2  2  (e2  -  1){1  +  d  (v  -  7)/a2}  +  o  (1); 

1=1  j=l  J  p 

A0Nl/2("»c(LL)  '  \)/a 

=  (NM)_1/2  2  2  (e2  -1)  +  d  M1/2N“1/2  2  (v  -  7)log  q2/a2  +  o  (1): 

1=1  j=l  v  1=1  1  i  v  p 

and 

AOn1/2^c(MML)  ’  ^c)/a 

=  (NH)~1/2  2  2  (e2  -  1)  +  dQM1/2N'1/2  2  (v  -  7)q2/a2  +  o  (1). 

1=1  j=l  u  1=1  1  i  v  p 

Simple  central  limit  theorem  calculations  yield  the  result. 

A 

Remark:  The  result  for  0^^  is  based  on  a  obtained  from  the  residuals  of  the 
final  fit  of  the  mean  response  function  as  for  the  log- 11  near! zed  method  and 
pseudo- likelihood,  so  that  Lemma  A. 3  holds.  The  modified  maximum  likelihood 
method  also  provides  a  Joint  estimate  of  a  along  with  the  estimate  of  0.  If 

A 

one  considers  this  estimator  in  place  of  a  in  Lemma  A. 3.  it  can  be  shown  that 
the  resulting  estimator  of  minimal  detectable  concentration  has  even  larger 
asymptotic  variance  than  that  in  Theorem  2. 

1 

Proof  of  Theorem  2  *•  By  (A. 5).  for  any  of  the  estimators  x  ,  since  q  =  £(x,p), 

c 

A  ^ 
we  have  for  some  A  that  £(x  ,p)  is  asymptotically  normal  with  mean  £(x  ,p)  and 

c  c 

2  A  M  1  /O  u 

variance  a  A/N.  Thus,  for  -t  between  x  and  x  ,  defining  WM  =  N  £  (x  ,/3 ) 

c  c  c  N  xy  c  ' 

(x  -  x  )  /  a.  we  have  WN  £(*»,/?)  /  t(x  ,/3)  is  asymptotically  normal  with 
zero  mean  and  variance  A.  where  *x(v.P)  is  the  derivative  of  the  first 
component  of  <(v,0).  It  thus  suffices  through  (A. 4)  to  prove  that 


£(t  ,P)/e(x  .P)  converges  in  probability  to  1.  But  this  follows  since  r>  / a 
c  c  c 

^(x_*0)^a  a_»  see  Lemma  A.l.  The  result  now  follows  from  Propostion  1. 

c  c 


Proof  of  Lemma  A .2-  By  a  series  of  Taylor  expansions  and  using  Lemma  A.l, 


N1/2{h(Sci)  -  h(0  ,?)}2/a2 

=  N1/2{h(T7*./3)  -  h(0 ,P))2/o2 

+  {d/drl  ^(O.P)}2N1/2(tjc  -  tj*)/ct  +  0  (1). 


1/0  ^  W 

Similar  calculations  noting  that  N  (P  -  P)  =  0  (a)  and  tj  -*  0  yield 

P  c 


N1/2h20(^  ,p)a2/(Ma2) 


=  {h20(7)^.p)/M)N1/2  +  (c/M)  h2B(0. P)N1/Z(a  -  a)/a2 

+  (2c/M)  h20(O.P)  {log  h(O.P)}  N1/2{0  -  0)  +  op(l). 


Combining  (1.5),  (3.1),  (A. 6)  and  (A. 7)  yields  (A. 5). 


Proof  of  Lemma  A. 3:  Define 


-12  I  !  2 


a  =  (NM)  ‘  2  2  [{Y  -  f(x  ,/3)}/fD(x..p)]z  =  (NM)'V2  2  2  ef, 

i=l  J=1  J  1  1  i=l  j=l  ij 


Then,  since  (NM )  (p  -  p)  =  Op(a). 


f*m\l/2/~2  ~2W  2 

(NM)  (a  -  an)/a 


N  M 


=(NM)  ^^{Yjj  -  f(x1./3))2/f20(xi.D) 


,  2  ,,20 


-  {Y  -f(x1.«)Vf",(x  WJo"-*  +  o  fl) 


20 


=  -2v(NM)1/2(0  -  0)  +  op(l). 


completing  the  proof. 
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Concentration  (x) 


Response  (Y) 

1.700.  1.660.  1.950,  2.070 


0.000 
0.075 
0.1025 
0.135 
0.185 
0.250 
0.400 
0.550 
0.750 
1.000 
1.375 
1.850 

I 

2.500 

3.250 

4.500 
6.000 

8.250 
•  11.250 

15.000 

20.250 

27.500 
37.000 


1.910.  2.270.  2.110.  2.390 
2.220  ,  2.250-,  3.260.  2.920 
2.800,  2.940,  2.380,  2.700 
2.780,  2.640,  2.710,  2.850 

3.540.  2.860.  3.150,  3.320 

3.910.  3.830,  4.880,  4.210 

4.540,  4.470.  4.790,  5.680 
6.060,  5.070.  5.000,  5.980 
5.840,  5.790,  6.100,  7.810 
7.310,  7.080,  7.060,  6.870 
9.880,  10.120.  9.220,  9.960 


11.040, 
13  10, 
16.070, 
17.340, 
18.980, 
21.666, 
23.206, 
23.922, 
25.748, 
24.441, 
29.580. 


10.460, 

15.470, 

14.670, 

16.850, 

19.850, 
21.218, 
22.239, 
24.871, 
25.874, 
25.874, 
26.698, 


10.880, 

14.210, 

14.780, 

16.740, 

18.750, 

19.790, 

22.436, 

23.815, 

24.907, 

25.748, 

26.536. 


11.650 

13.920 

15.210 

16.870 
18.510 
22.669 
22.597 

24.871 
24.871 
27.270 
27.181 


50.000 


*PL 

*LL 

*MML 

xc(LS) 

xc(PL) 

xc(LL) 

xc  (MML) 

Full 

Duplicates 

0.4750 

0.4757 

0.4500 

0.1554 

0.0790 

0.0793 

0.0822 

1  k  2 

0.7000 

0.9404 

0.7500 

0.2230 

0.0728 

0.0476 

0.0659 

2  &  3 

0.3500 

0.1950 

0.2500 

0.2385 

0.1555 

.  0.1870 

0.1739 

3  &  4 

0.5750 

0.6940 

0.7000 

0.2513 

0.1324 

0.1112 

0.1104 

1  &4 

0.5500 

0.5931 

0.5000 

0.1593 

0.0612 

0.0601 

0.0695 

1  k  3 

0.4500 

0.4233 

0.4750 

0.1859 

0.0938 

0.0981 

0.0909 

Mean 

0.5250 

0.5692 

0.5350 

0.2116 

0.1031 

0.1008 

0.1021 

SD 

0.1183 

0.2511 

0.1997 

0.0763 

0.0357 

0.0491 

0.0439 

Note:  Means  and  SDs  are  based  only  on  the  five  reduced  combinations  of  the  data  with  duplicates. 
Estimate  of  minimum  detectable  concentration  based  on  unweighted  least  squares  is  denoted  xc(LS). 


t 


